Non-abelian gauge fields and quadratic band touchings in molecular graphene 
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Dirac fermions in graphene can be subjected to non-abelian gauge fields by implementing certain 
modulations of the carbon site potentials. Artificial graphene, engineered with a lattice of CO 
molecules on top of the surface of Cu, offers an ideal arena to study their effects. In this work, we 
show by symmetry arguments how the underlying CO lattice must be deformed to obtain these gauge 
fields, and estimate their strength. We also discuss the fundamental differences between abelian and 
non-abelian gauge fields from the Dirac electrons point of view, and show how a constant (non- 
abelian) magnetic field gives rise to either a Landau level spectrum or a quadratic band touching, 
depending on the gauge field that realizes it (a known feature of non-abelian gauge fields known 
as the Wu-Yang ambiguity). We finally present the characteristic signatures of these effects in 
the site-resolved density of states that can be directly measured in the current molecular graphene 
experiment, and discuss prospects to realize the interaction induced broken symmetry states of a 
quadratic touching in this system. 



I. INTRODUCTION 

Condensed matter systems that host Dirac fermions 
as their electronic excitations have drawn a lot of atten- 
tion in recent years as they have become more and more 
experimentally accesible and controlable, with graphene^ 
and topological insulators^ being the most prominent ex- 
amples of such materials. 

A remarkable feature of Dirac fermions realized in 
graphene's honeycomb lattice in particular is that one 
can further manipulate them externally by inducing con- 
trolled strains in the sample, which couple to them 
as an effective gauge potential 3 . This idea of strain 
engineering^ has led to many interesting predictions^, 
and is most spectacularly illustrated by the Landau level 
spectrum recently observed^ in scanning tunneling mi- 
croscopy (STM). This system proved to be very versatile, 
and in search for even better tunability several proposals 
were conceived to make artificial versions of itPESl. In 
a recent experimental breakthrough, a realization of this 
type of systems, termed molecular graphene 11 , was built, 
which allows for almost complete control of the electronic 
degrees of freedom within it. In this system, a triangu- 
lar lattice of CO molecules is assembled in the surface 
of bulk Cu, confining the surface electrons to move in 
an effective hexagonal potential. In this way, effective 
Dirac fermions emerge at the K points of the superlat- 
tice potential, which can then be probed directly with an 
STM. 

This system thus offers wide tunability to modify the 
electronic structure of the surface states by distorting the 
CO lattice in any desired way, or by adding new atoms to 
the existing structure. Indeed, several remarkable phe- 
nomena have already been demonstrated^ beyond the 
strain induced Landau levels, such as the opening of a 
gap by means of a Kekule distortion or the creation of 
an n-p-n junction. Other interesting proposals such as 
the observation of fractional charge in a vortex*^^ or 
the synthesis of a quantum spin Hall phase^ should also 



be experimentally accesible. 

As noted in ref. [15j artificial graphene should be also 
ideal to explore the more recent prediction that a full 
SU (2) non-abelian gauge field is in fact realizable in this 
system, and the strain induced one is just one component 
of it. Non-abelian gauge fields (of singular nature) were 
known to emerge in graphene due to disclinations in the 
lattice^, but they can also be generated in a smooth 
fashion by modulating the on-site potential of the carbon 
atoms in a certain way. As we will discuss, the effects of 
non-abelian gauge fields can be very different from their 
abelian counterparts, and it is the purpose of this work to 
discuss how to adapt the molecular graphene experiment 
to probe these differences. In particular, we will show 
that a quadratic band touching can be generated with 
these fields, allowing a controlled simulation of this band 
structure which is prone to many-body instabilities^^. 

In general, effective external gauge fields acting on a 
fermion system may have non-abelian structure when the 
fermions have internal degrees of freedom, and the gauge 
field is a matrix acting on this degree of freedom 
whose components need not commute. A typical con- 
densed matter example is spin and the spin-orbit interac- 
tion, which can be modeled as an SU(2) gauge fiel d 20 * 21 *, 
but there are many more examples^H 2 ^. A more recent 
one is bilayer graphene^, where the two components 
of the SU(2) doublet correspond to the wave functions 
in the two layers, and the interlayer interaction plays 
the role of the gauge field. In the case of monolayer 
graphene, the SU(2) doublet is made with the valley de- 
gree of freedom^. 

The non-abelian field strength is defined in terms of the 
covariant derivative Di = di — iAi as Fij = [D^Dj\ = 
diAj — djAi — i[Ai,Aj), which in two dimensions gives 
rise to a non-abelian magnetic field of the form 

B a = dxA a + e aM A^ x A 7 , (1) 

where A ah = A a A% b with A% b the generators of SU(2), re- 
peated indices are summed, and a = x,y,z (the indices 
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TABLE I. Classification of basis matrices in the low energy 
theory around the K,K' points in graphene according to the 
representations of the symmetry group Cq v , and their explicit 
realization in the basis (^ak^bk^ak'^bk 1 ) (see ref. [30] 
for details). 



ab will be implicit from now on) . The last term in this ex- 
pression arises because of the non-commutativity of the 
field components and makes non-abelian gauge fields fun- 
damentally different from their abelian counterparts. In 
particular, it is responsible for a tricky feature of these 
gauge fields known as the Wu-Yang ambiguity^ the fact 
that one may have physically distinct gauge fields (i.e. 
not gauge equivalent) with the same magnetic field. In- 
deed, consider these two simple examples 28 . The first 
(type I) is = B/2(-y,x), A« = A^ = 0, which 
we recognize as the analog of the symmetric gauge for 
constant (abelian) magnetic field B, in this case in the 
z direction. The second (type II) is A^ = ^5/2(1,0), 
= y/W/2(0, 1) and A^ = 0, it also gives constant 
field Bo due to the second term in eq. ([!]), and it is not 
gauge related to type I. The magnetic field alone is there- 
fore not enough to distinguish these two cases^, but we 
will see that the spectrum obtained for each case is very 
different, and this is the physics that, as we will show, 
can be probed directly in the molecular graphene exper- 
iment. 



II. SYMMETRY ANALYSIS AND 
MICROSCOPIC CALCULATION 

To realize an SU(2) gauge field in graphene, what 
we need is to externally apply certain on-site potential 
patterns^ to the carbon atoms of the honeycomb lat- 
tice. It is not a priori clear, however, how this may be 
achieved in a molecular graphene experiment, where the 
"effective honeycomb lattice" is engineered with the po- 
tential landscape induced by a triangular array of CO 
molecules. In terms of an effective tight binding model, 
it is natural to think that small distortions of this tri- 
angular lattice will produce potential changes in the ef- 
fective carbon sites, but what distortions will give rise 
to the correct potentials? And more importantly, since 
these distortions may induce changes in the effective hop- 



ping as welffH, is it possible to modulate only the on-site 
potential? 

To answer these questions, a symmetry approach to 
the problem appears better suited. The way that ex- 
ternal perturbations couple to the low energy degrees of 
freedom around a high-symmetry point of the Brillouin 
Zone can be determined just by symmetry arguments. 
This approach has been fruitfully employed in graphene 
to discuss the coupling of phonons, strains, or electro- 
magnetic fleld d 30 " 33 ! and we now show how it can be used 
to see the emergence of non-abelian gauge fields from 
small CO displacements the molecular graphene. 

In the half- filled honeycomb lattice, electrons close 
to the Fermi surface live near the K and K' 
points, and are described by an effective spinor 
(ipAK 5 ip bk i i^AK' i *0 BK' ) •> where A/B denotes the sub- 
lattice degree of freedom. The effective Hamiltonian is 
conventionally written in the basis of the Pauli matrices 
o~i (g) Tj, where acts on the sublattice and on the 
valley degrees of freedom, and i = x,y, z (the identity 
in both sets is understood to be included as part of the 
basis). 

To exploit the fact that the Hamiltonian must be a 
scalar under the symmetry group Cq v of the honeycomb 
lattice, one can relabel these basis matrices in terms of a 
new symmetry adapted set E^ and A$ with the Pauli ma- 
trix algebra and well-defined transformation properties 
under this group (technically, the group is C$ v because 
the unit cell has been tripled to consider K and K' at 
the same time. We will refer to the labels under Cq v for 
simplicity, see ref. |30| for details) . The relation of these 
matrices to the original ones and the representations ac- 
cording to which they transform are reproduced in table 

m 

In this basis, the low energy Hamiltonian is simply 
written as (vp = 1) 

H = E • fc, (2) 

and in this form it is simple to see that the matrices A^ 
commute with the Hamiltonian and generate an SU(2) 
symmetry, which corresponds to rotations in the valley 
degree of freedom. 

A gauge field is by definition a field that couples mini- 
mally in the form ki —> ki + Ai, and in analogy with the 
usual electromagnetic field that couples as Hum — E- A, 
one may introduce an SU(2) gauge field that couples as 

H su(2) = S • (a x AW + A y AM + -K Z A^) , (3) 

which is a coupling allowed by symmetry if the gauge 
fields A a have their origin in a microscopic perturbation 
with the same symmetry as the matrix that accompanies 
them. 

The power of the symmetry analysis is thus that one 
can now say what type of perturbations correspond to 
each term only by inspection of table [I] Perturbations 
in the first column have the periodicity of the unit cell, 




FIG. 1. (Color online) The four possible CO displacements 
with symmetries E\ and E2. The CO molecules are repre- 
sented in red, and the effective honeycomb lattice is shown in 
black. The unit cell is shaded in gray, but more hexagons are 
shown to make the symmetry of the modes apparent. The 
on-site potentials that match the symmetry labels are also 
shown in the effective carbon sites. Note that the pref actors 
only refer to the potentials, not to the displacements. One 
may think of these displacements as the K-point phonons of 
the triangular CO lattice. 



while those in the second column have the periodicity of 
a tripled unit cell (because of intervalley mixing) . More- 
over, within nearest neighbour tight binding (TB), those 
perturbations diagonal in sublattice (ex ao or a z ) corre- 
spond to potential modulations, while those off diago- 
nal correspond to hopping modulations. With this crite- 
rion, the gauge field is readily identified as the usual 
strain-induced gauge field. The gauge field components 
A( x ) and A( y ) correspond, respectively, to the valley mix- 
ing Ei and E2 potential perturbations defined in ref. [30] 
(which are labeled G' under Cq v ). Their corresponding 
potentials are depicted in fig. [I] in the effective carbon 
sites. In real graphene this type of potential perturbation 
is the one induced by phonons like the LO/LA phonon 
at the K point 34 , or the ZO/ZA phonon at the K point 
in the presence of a perpendicular electric field, and it 
is known that it can also be produced by a particular 
substrate 35 . 

For molecular graphene, this analysis immediately al- 
lows to find the CO displacements that will induce these 
potential modulations. These should be displacements 
with a tripled unit cell and the appropriate symmetry 
labels, and in fact may be simply interpreted as the Ei 
and E2 phonons of the triangular CO lattice at the K 




FIG. 2. (Color online) With the same conventions that fig. 1, 
combinations of CO displacements that produce a quadratic 
band touching. Again, note that the prefactors only refer to 
the potentials, not to the displacements. 



point. These displacements are, for three consecutive 
CO molecules 

rco,E lx = {(1,0), (-|,0), (-|,0)}, 
fco, Ely = {(0,1), (0,-1), (0,-|)}, 

fco,E 2x = {(0,0), (O,-^ 3 ), (O,^ 3 )}, 

W 2 , = {(o,o), (^,0), (-^,0)}, 

and are also shown in fig [I] Indeed, within a TB model 
one can parametrize the change in on-site potential with 
displacement as 

V i = V f J2 A ^co-S ij , (4) 

3 

for a carbon site i with j CO neighbours at equilibrium 
distances 5{j from it, and verify that the potentials shown 
in fig. [I] are given by eq. Q. The constant V = dV/da 
parametrizes the change in on-site potential with dis- 
tance, and may be estimated by realizing that this phys- 
ical mechanism is responsible for the scalar potential <fi 
in the continuum Dirac equation. A comparison with 
the p-n junction experiment yields V « 22 meV/A (see 
appendix), which is very similar to dt/da = fit /a « 20 
meV/A. This is also consistent with the fact that in real 
graphene the analog of V for carbon displacements^ is 
of the same order as dt/da. 
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Finally, the symmetry analysis also reveals that close 
to the Dirac point, the desired CO displacements do not 
introduce any other change in the effective theory other 
than the A^ X \A^ gauge fields. In particular, while 
these displacements may induce nearest neighbour hop- 
ping changes, these cannot appear in the low energy the- 
ory because there are no intervalley matrices in the E\ or 
E2 representations that are sublattice off-diagonal. This 
hopping changes thus have no effect in the low energy 
properties and we will not consider them in what fol- 
lows. Changes in the next nearest neighbour hopping t' 
due this displacements are small and need not be consid- 
ered. 

To obtain the gauge field from a microscopic calcula- 
tion, one may substitute eq. Q in the effective tight 
binding model 

H = ~ l E c l°j - l ' E c l°i + E M*. (5) 

(hj) «M» i 

The potential modulation (depicted in fig. [I]) that gives 
rise to the non-abelian gauge fields iJ^ 

3 r 

V(x) = -V' u E2y cosKx + u Elx sinKx (6) 



l Gx (u Ely 



cos Kx + ue 2x sin Kx 



with K = (4tt/3V3, 0) a vector joining the two dirac 
points and G = (0, — 47r/3) a reciprocal lattice vector. 
To project this perturbations into the Dirac points one 
performs the sum 

H = Y.Md = Y, V (^4 + V(x + 5 1 )c^c?, (7) 



with x = ndi + ma2 the lattice positions, Si — a(0, 1) a 
nearest neighbour vector, and 



A _ iKx A 1 —iKx A 
c x — e C K * e C K' •> 



JKx B 
C K 



iKx^B 

c K , 



(8) 
(9) 



This sum gives exactly the matrices dictated by symme- 
try 

3 

H= ^V , (-T2UE lx -\-r 1 a 3 UE ly -r 2 (TsUE 2x -\-r 1 UE 2y ), (10) 

so that the final formula relating the effective gauge fields 
to CO displacements is 



AW = -V'(-u Ely ,u El J, 
= \V{u E ^,u E2y ). 



III. PHYSICAL EFFECTS 



(11) 
(12) 



As described in the introduction, we now consider two 
^auge field configurations that are not related by a gauge 




FIG. 3. (Color online) Total density of states for any type II 
gauge field of strength u = 0.5 A (red line) and u = lA (black 
line), for t = 90 meV and t' — 0. The unperturbed LDOS 
is shown as a dashed blue line for comparison. Inset: Band 
structure of the system for u — lA. Note the similarity with 
bilayer graphene. 



transformation, but whose magnetic field is the same, 
and consider how they should be seen in a local den- 
sity of states (LDOS) measurement. Consider the type 
I gauge field, with a magnetic field pointing in a general 
direction b a in 517(2) space, Af = b a B/2(y, -x). When 
b a = (0, 0, 1) we have the usual strain induced gauge 
field. The case b a = (0,1,0) was discussed in ref. [T5l 
In general, by a constant SU(2) rotation of the Hamil- 
tonian, it is not difficult to see that for any b a the spec- 
trum is still given by Landau levels E n — \JlBn. The 
only difference appears in the wavefunctions, because the 
sublattice polarization turns out to be given by the pro- 
jection of b a onto the z axis. For strain induced fields 
it is maximum, but for potential induced ones the den- 
sity of states is in fact constant across the unit cell. One 
can estimate the magnetic field induced in the molecular 
graphene experiment with these gauge fields as follows. 
Take u Elx = u max /L * x, with L the radius of the (ap- 
proximately circular) sample and u max the maximum dis- 
placement (at x = L). The magnetic field is (recovering 
all units) 



H/e W 



L 



(13) 



With h/e = 6.5 • 10 4 TA and taking u max = 0.1a and 
y/3a/L ^1/10 and a Fermi velocity 11 hvp ~ 1.5eVA we 
obtain B w 3.75T, which is not very large compared to 
the strain induced one that is typically achieved. 

The type II gauge field has better prospects to be ex- 
perimentally accessible. Keeping A^ = 0, there are in 
fact four possible choices of constant gauge fields that 
give constant magnetic field, given by 

A 1 = y/B/2{l, 0) A 2 = VW2(0, ±1), (14) 
which, by eq. (|12|), is produced with the displacements 
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E/t E/t 



FIG. 4. (color online) LDOS as a function of energy for two 
different gauge fields. Top left: E\ x — E 2x , with strength u=l 
A, t = 90 meV and t' — 0. Top right: E\ y — E 2y , same 
parameters. Bottom plots are the same but with t 7 = 0.18t 
and a Lorentzian broadening of E = 0.2t. The insets show 
the corresponding on-site potential and the color code for the 
different sites within the unit cell. Note that the missing lines 
in the plots overlap with the line shown that has the same on- 
site potential. 



Ely ± E 2y , and 

A 1 = y#72(0, 1) A 2 = y/B/2(±l, 0), (15) 

which is produced with the displacements E\ x ± E 2x . 
These displacements and their on-site potentials are de- 
picted in fig. [2] The magnetic field is given by B = 
9/8(V'u/vf) 2 with u = ue xi = ±UE 2i representing the 
modulus of the displacements in fig. |2j It is interesting 
to note that the estimate for B in this case for u = 0.1a 
is B « 16T. 

The Dirac Hamiltonian in the presence of these gauge 
fields is formally analogous to that of bilayer graphene 
(for a single valley), with the role of the layer played 
by the valley here^, and an effective interlayer coupling 
7 = \[2B. The spectrum of these Hamiltonians is well 
known to be a quadratic band touching, with two extra 
parabolic bands at higher energies. 

Considering first the case t' = 0, the density of states 
(DOS) of this system is finite at the touching point Ed = 
0, and has a jump at ±vfV%B = 3/2V'u, as depicted 
in fig. [3J Considering a displacement u = 0.1a = lA, 
the kink in the LDOS should appear at ±30 meV, which 
should be easily observable. The precise location of this 
jump should serve as an independent estimate of the pa- 
rameter V . The main effect of a finite t' is to shift Ed 
to a higher value, as we will see below. 

Moreover, this type of gauge field shows more compli- 
cated local density of states across the enlarged unit cell. 
In fig. [4] we show the LDOS for the cases E\ x — E 2x and 
Ely — E 2y . The other two combinations are obtained by 



mirror symmetry. For t' = 0, we observe different lo- 
cal gaps for different sites, and finite LDOS at E = 0. 
For more faithful comparison with the experiment, we 
have also plotted the LDOS for t' = 0.18t, and with 
a Lorentzian broadening of E = 0.2t. We observe the 
main effect of a shift in Ed, as well as some electron- 
hole assymmetry, but the main features that characterize 
the non-abelian gauge field remain. The identification of 
these features in an STM measurement would represent 
a demonstration of the presence of the type II constant 
non-abelian gauge field. 



IV. DISCUSSION 

In this work we have shown, by means of a sym- 
metry analysis, how non-abelian gauge fields may be 
implemented in molecular graphene, and what their 
experimental signatures should be in the LDOS. For 
type I gauge fields of constant magnetic field, we have 
shown that because of the different microscopic origin of 
gauge fields A^ (hopping change) and A^ 1,2 ^ (potential 
change), the magnetic field that one gets in the second 
case is relatively smaller. While this may make the Lan- 
dau level spectrum more difficult to observe, the presence 
of this type of field could also be readily detected, for ex- 
ample, in a quantum interference experiment in the weak 
field limiP. 

We have also shown that type II constant non- 
abelian gauge fields generate a quadratic band touch- 
ing analogous to bilayer graphene. Because of the en- 
hanced DOS at the Fermi level, the electron-electron 
interaction is known to drive this system to a broken 
symmetry state whose precise characteristics are still 
controversia l 17 " 19 ^ . In the current molecular graphene 
experiment, the Coulomb interaction is screened by the 
metallic bulk, leaving residual Hubbard interactions es- 
timated to be U ~ 0M ~ 50 meV (see Supplementary 
Material of ref. ITT]) . While an ideal quadratic touch- 
ing is unstable to infinitesimal short range interactions, 
the current broadening due to bulk tunneling (~ 0.2t) is 
perhaps too large and may challenge the observation of 
the interaction induced transition. Both bulk tunneling 
and screening could be reduced by performing future ex- 
periments in bulk insulators with metallic surfaces (such 
as the recently discovered topological insulators^) which 
may eventually allow to study the fate of the many body 
state with a tunable analog of the interlayer hopping. 
Incidentally, it is also interesting to note that this in- 
stability can be interpreted as a non-abelian magnetic 
catalysis, where an infinitesimal field drives chiral sym- 
metry breaking 37 . Furthermore, the controlled simula- 
tion of these non-abelian gauge fields, when made posi- 
tion dependent, may be used to study the generation of 
zero-energy flat bands, as those observed in the twisted 
bilayer systenP|, or the physics of topological defects in 
the gauge nelcPn 

In summary, the molecular graphene experiment has 
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great potential to observe many interesting phenomena 
related to non-abelian gauge fields with an unprecedent 
tunability, and which, as we have shown, should be real- 
izable in the current experimental samples. 
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VI. APPENDIX 

A. Estimate of V 

The parameter V describes the change of on-site po- 
tential due to the displacements of neighbouring CO 
molecules. As such, it is featured both in the non abelian 
gauge fields (which come from "optical" displacements) 
and in the strain-induced scalar potential <j) (which comes 
from "acoustical" displacements). The scalar potential (j) 
also has a contribution from NNN hopping change dt' /da 
but it is much smaller and will be neglected. 

To see this, consider an isotropic expansion of the CO 
lattice. For every carbon site z, the induced potential is 
given by eq. Q. Because displacement is smooth we 
may write 



s L djri Xj co 



(16) 



and plugging directly into the TB Hamiltonian eq. ([5]), 
we obtain that <j) = 3V'a/2(u xx + u yy ). Now consider the 
p-n-p juntion experiment in ref. [Til The middle region is 
strained from d = 17.8 to d = 20.4 so u xx = u yy = 0.14. 
The change in scalar potential A<p is 95 meV, so we obtain 
(d = y/3a) 



v 



95meV 



O.Uy/S 17.8A 



22 meV/A. 



(17) 



A different estimate can be obtained from the nearly free 
electron model considered in ref [TT1 (supp. mat.), where 



the scalar potential is 
8tt 2 



which gives V = 24 me v/A 

B. Symmetry of hopping perturbations 

There are 9 independent hoppings in the tripled unit 
cell, which can be decomposed into combinations that 
have well defined transformation properties under the 
symmetries of the lattice. The 9 combinations and their 
symmetry labels are shown in fig. [5] In the first row, one 
may identify the constant hopping (Al), the E2 pattern 
that gives rise to the usual gauge field, and the Kekul 
distortions (any of the three domains can be obtained 
from these). 

The four combinations in the second row are E± and E2 
(and form the representation C when the enlarged group 
Cq v is considered) . These hopping patterns are produced 
by the same CO displacements that give the non-abelian 
gauge fields through charge modulation. In the main 
part of the text we claimed that these hopping distor- 
tions cannot couple to the low energy theory around the 
K point. The reason is simply that there is no valley off- 
diagonal Ei or E2 matrix in the low energy theory whose 
microscopic origin is a hopping change. This can be seen 
directly by inspection of table [T[ where the valley mixing 
Ei and E2 matrices are all diagonal in sublattice. The 
hopping modulations will only appear in the low-energy 
theory if terms with higher order in momentum are con- 
sidered. If one is interested in the whole band structure 
and not just low energies, these distortions in the hop- 
ping should be included by changing the NN hopping in 
the usual manner. 




E 2x 



Ely 




E2.X 



E2 y 



FIG. 5. (Color online) The 9 independent hopping patterns 
and their symmetry labels. Blue means positive and red nega- 
tive, and black lines represent no change in the hopping. Hop- 
ping modulations of the corresponding symmetry may also be 
induced by these displacements (red is negative hopping and 
blue is positive. However, these particular patterns have no ef- 
fect in the low energy theory: they do not affect the quadratic 
touching and the LDOS predictions around E — Ed- 
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